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ABSTRACT 

Recently, QED bremsstrahlung in 5-meson decays into pair of scalars (Tt's and/or ^"s) is 
of interest. If experimental acceptance must be taken into account, PHOTOS Monte Carlo is 
often used in experimental simulations. We will use scalar QED to benchmark PHOTOS, even 
though this theory is of limited use for complex objects. We present the analytical form of 
the kernel used in the older versions of PHOTOS, and the new, exact (scalar QED) one. Matrix 
element and phase-space Jacobians are separated in the final weight and future extensions based 
on measurable electromagnetic form-factors are thus possible. 

The massive phase-space is controlled in the the program with no approximations. Thanks to 
the iterative solution all leading and next to leading logarithmic terms are properly reproduced 
by the Monte Carlo simulation. Simultaneously, full differential distributions over complete 
multiple body phase-space are provided. An agreement of better than 0.01% with independent 
calculations of scalar QED is demonstrated. 
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1. INTRODUCTION 



In the analysis of data from high-energy physics experiments, one tries to resolve the "ex- 
periment = theory" equation. This non-trivial task requires that a lot of different effects be 
considered simultaneously. From the experimental side, these are mainly detector acceptance 
and cuts, which are dictated by the construction and physical properties of the detector. The 
shapes of distributions may be distorted by, say, misidentification and residual background con- 
tamination. These effects need to be discriminated in an appropriate and well-controlled way. 
From the theoretical side, all effects of known physics have to be included in predictions as 
well. Only then experimental data and theoretical predictions can be confronted to determine 
numerical values of coupling constants or effects of new physics (to be discovered). 

A well-defined class of theoretical effects consists of QED radiative corrections. PHOTOS 
is a universal Monte Carlo algorithm that simulates the effects of QED radiative corrections in 
decays of particles and resonances. It is a project with a rather long history: the first version 
was released in 1991 [1], followed by version 2.0 [2] (double emission and threshold terms for 
fermions). The package is in wide use [3]; it was applied as a precision simulation tool for the 
W mass measurement at the Tevatron [4] and LEP [5,6], and for CKM matrix measurements in 
decays ofK and B resonances (NA48 [7], KTeV [8] , Belle [9], BaBar [10] and Fermilab [11]). 
Discussion of the different components of systematic errors in PHOTOS is thus of interest. 

Throughout the years the core algorithm for the generation of 0(a) corrections did not 
change much, however, its precision, applicability to various processes, and numerical stability 
improved significantly. New functionalities, such as multiple photon radiation and interference 
effects for all possible decays were introduced [12, 13]. Recently, the complete first order ma- 
trix element was introduced into PHOTOS for Z decays and complete NLO0 multiple photon 
predictions for that channel were demonstrated to work well [14]. 

Increasing interest in the algorithm expressed by experimental collaborations (including 
future LHC experiments and precise measurements for B decays) was a motivation to perform 
a more detailed study of the potential and precision of the PHOTOS algorithm. This paper is 
devoted to the decay of 5-mesons into a pair of scalars. It is a continuation of the previous paper 
[14] devoted to Z decays. We concentrate our attention on exact phase-space parametrization as 
used in PHOTOS, and on the explicit separation of the final weight into parts responsible for the 
following: (i) mass dependent phase-space Jacobians, (ii) matrix elements and (iii) pre-sampler 
for peaks. 

Simplifications introduced in the matrix element normally used in these scalar B meson 
decays are removed, the exact kernel of first order scalar QED calculation is installed. Such 
improvement opens the way to include data dependent form-factors into matrix elements of 
PHOTOS, and physically better results than of scalar QED alone. 

'in the paper, we will use abbreviations NLO, NLL, NNLO, NNLL to denote next to leading order, next to 
leading logarithm, next-next to leading order, next-next to leading logarithm corrections with respect to leading 
order (that is without QED at all). Meaning for such abbreviations can not be defined restrictively on the basis of 
approximations used in phase-space. The properties of the matrix elements need to be specified as well. Nonethe- 
less, we will use the abbreviations later in the paper to denote such approximations in the phase-space Jacobians 
which do not prevent the appropriate precision to hold, if proper choice of matrix elements is made as well. 
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Our study of the PHOTOS algorithm can be understood as another step in the on-going 
effort to find practical solutions of the improved expansions. The solution can be seen as a 
rearrangement of the QED perturbation expansion, but this time for the interaction of charged 
scalars with photons and in the case where ultrarelativistic approximations are not valid. 

To test PHOTOS we have used predictions of the SANC [15] Monte Carlo algorithm. SANC 
is able to calculate the exact first order scalar QED matrix elements for decays of 5-mesons 
into scalars, and covers the full phase space of decay products without any approximations. 
Events provided by SANC MC are unweighted. SANC is a network client-server system for 
the semi-automatic calculation of Electroweak, QCD and QED radiative corrections at a one- 
loop precision level for various processes(-decays) of elementary particle interactions. 

The paper is organized as follows. Section 2 is devoted to the description of our results 
obtained from scalar QED, which will be used later in tests and in construction of kernel for 
single photon emission. In Section 3 the main properties used in the design of PHOTOS are 
presented. In particular, the construction of the weight (NLO level) necessary to introduce the 
complete first order matrix element is explained in all detail. The phase space parametrization 
used in the iterative solution of PHOTOS is given. Details are collected in the Appendix. On the 
other hand, quite essential for complete NLL, aspects of the construction will still remain poorly 
documented. That is, how the parts of the single photon emission matrix element are used in the 
(iterated and thus extended to the mutiphoton emissions) kernel. Section 4 is devoted to results 
of numerical tests performed at fixed, first order of QED. Some numerical results, obtained 
with multi-photon version of the program will be collected there for reference. Finally, section 
5 summarizes the paper. 



2. Scalar QED and B decays 



In this section we give the formulae needed in the construction of the Monte Carlo routine 
of two-particle 5-meson decays and the analytical results of the decay rates at O(oc) in scalar- 
QED when the masses of the decay products are neglected. The one-loop QED corrections to 
the width of the decays B°~ — > H~H+ °, where H l 2 denotes scalar(pseudo-scalar) particles, 
can be represented as a sum of the Born contribution and the contributions due to virtual loop 
diagrams and soft and hard photon emissions. Both virtual and soft contributions factorize to 
the Born one. 
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Here jr Born is the tree level differential decay width, h virt represents the virtual corrections, 
§«>/f denote j ne so ft photon contribution and jr Hard is the hard photon contribution. The Born 
level distribution in the rest frame of the decaying meson can be written as 
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where M is the mass of decaying particle, k\^ denote the momenta of decay products, A Born 
stands for the corresponding tree level amplitude and dLipsiiP — > k\,k2) is the two body dif- 
ferential phase space. For the latter we choose the following parametrization 
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where angles 0i and (|)i define the orientation of momentum k l in the rest frame of the 5-meson. 

As usual we define X^{a,b, c) = \J a 2 + b 2 + c 2 — lab — lac — Ibc. In the case of neutral B- 
meson decay channels B° — > H~H 2 + , the scalar-QED calculations for the virtual and soft factors 
in formula (OQ) gives 
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where m, 2 are the final meson masses, an auxiliary small parameter a)(<^ M/2) separates the 
soft and hard photon contributions, and fi uv denotes the ultraviolet scale. An auxiliary photon 
mass my is used as a regulator of the infrared divergences. The ultraviolet singularities are 
regularized by means of the dimensional regularization. We renormalize the wave functions of 
the external scalar fields in the on-shell scheme, and the point-like weak coupling in the MS 

^ 2 -i 2 ) andLi 2 (z) = - T — In 1 1 
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The hard photon distribution dT 
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in scalar-QED can be expressed as follows 
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Here, are the charges of final mesons, and ky and are the photon momentum and po- 
larization vector respectively. The three body differential phase space of the decay products, 
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dLipsT,(P — > k\,k2,ky), is parametrized, in a rather standard way (see e.g. [16]), as follows: 
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where the angles 9f and § R define the orientation of momentum k\ in the rest frame of (ki + 
ki\, the photon energy Ey and the angles 9y and <|>y, that define the orientation of the photon 
momentum, are given in the rest frame of the decaying particle. These parameters vary in 
the limits: < Y , 0f < ji , < 4> Y , 4>f < 2n and CO < Ey < (M 2 - (mi + m 2 ) 2 ) /2M. The 
formulae © - © are used in the construction of the Monte Carlo simulator of the decays under 
consideration. An analytical integration in © over the phase space variables © can be done 
easily. Below we give the result of integration in the massless limit of the final mesons (i.e. 
m ll m 2 = m — > 0) because of its simplicity: 



^Hard 



pBorn^ 
71 



m 



1 — In — ~ In 

M 2 / M 2 



4eo 2 „, m 2 7t 2 

21n— T h4 

M 2 3 



(8) 



The virtual correction depends on the ultraviolet scale /j uv , which should cancel in the total 
decay width because of the scale dependence of the point-like weak coupling. The infrared 
divergence cancels out in the sum of virtual and soft contributions, as it must. The total decay 
width, which is the sum of the contributions ©, ©, and ©, is also free of go and of the final 
meson mass singularity in accordance with the KLN theorem [17]- [18]: 
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The same calculations can be done for the charged 5-meson decay channels B 
In this case for various contributions in formula (OQ) we obtain: 
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- Soft photon contribution 
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Again after integration over the phase space variables, the massless limit of the final mesons 
(i.e. m 1 ,m 2 = m — > 0) was used in the last formula. Finally, summing contributions (flOl) . (fTT|) 
and (fl"2~l) . we obtain the following expression for the total decay width: 
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We have checked that the factors 8 soft and 8 virt in formula (QQ), for both charged and neutral 
5-meson decays, provide the same numerical results as the corresponding expressions in [19]. 

To be assured of the accuracy of SANC Monte Carlo integration (which is a by-product of 
MC simulation), we compared Monte Carlo results with the analytical calculations of the total 
decay rate. The result of this comparison is shown in Table lH. The agreement is thus better 
than 10~ 4 in this test, where mass effects were included. 



Channel 


T^ al , lO^MeV 


T^ al , 10~ 3 MeV 


B~ -»• 7C-7C 

B -> K~K° 


0.373629 
0.367586 


0.3736(4) 
0.3675(9) 


B U -> 7C - 7l + 

5° -> K~K+ 


0.377392 
0.371414 


0.3773(8) 
0.3714(2) 



Table 1 : Comparison of the B-meson total decay widths produced by analytical calculations 
( second column ) and by Monte Carlo ( third column ). The last significant digit of the Monte 
Carlo results is given in bracket. 



-Please note that these numbers are for the purpose of our test only, the overall B—H — H coupling constants 
do not match the experimental data. 
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3. Exact phase-space and matrix element. 



To start any discussion of the implementation of complete first order QED radiative correc- 
tions in B decay, one has to specify the parametrization of the complete phase-space slots of the 
fixed final state multiplicity. 

Let us start with the explicit expression for the parametrization of an n + 1 body phase- 
space in decay of the object of four-momentum P (P 2 = M 2 ), as used in PHOTOS Monte 
Carlo. As our aim is to define iterative relations, let us denote the four momenta of the first n 
decay products as fc, (i = 1 , n) and the last n + l decay product as k n+ \ . In our case the n + 1-th 
particle will always be the real and massless photoro In the later steps of our construction the 
masslessnes of photons and properties of QED matrix elements will be used. 

In the following, notation from refs. [20, 21] will be used. We will not rely on any particular 
results of these papers. We only point to other, similar options for the exact n-body phase-space 
parametrizations, which are also in use. 

The Lorentz invariant phase-space is defined as follows: 
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where extra integration variables, four vector p (compensated with 8 4 (p — YJ[ h) ) is introduced. 
If further, Mi... M (compensated with d(p 2 — M 2 n ) ) is introduced, the element of the phase space 
takes the form: 

dM 2 

dLips n + 1 (P) = . L - n dLips 2 (P -> p k n+i ) x dLips n {p^k\...k n ) 
(2tc) 



dM 2 
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x dLips n (p — > k\ . . .k n ). 

(15) 



The part of the phase space Jacobian corresponding to integration over the direction and en- 
ergy of the last particle (or equivalently invariant mass Mi... n of the remaining system of \...n 
particles) is explicitly given. We will use later in the formulas m 2 = k 2 , and analogously A/, . 
defining invariant masses of k[ ...k n systems. The integration over the angles and (]) is defined 
in the P rest-frame. The integration over the invariant mass, M\ is limited by phase space 



3 However the construction does not rely on a photon to be massless. In principle it can be applied to define 
other phase space relations, for example the emission of an extra massive pion or emission of a pair of heavy 
particles. 
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boundaries. Anybody familiar with the phase-space parametrization as used in FOWL [16], 
TAUOLA [21], or many other programs will find the above explanation quite standard^. 

The question of choice of axes with respect to which angles are defined, and order in kine- 
matical construction, is less trivial. The choice for the particular option stems from necessity 
to presample collinear singularities. It is rather well known that the choice of the reference 
directions for the parametrization of the unit sphere is free, and can be used to advantage. We 
will use related, but somewhat different freedom of choice. Instead of variables <j) defining 
orientation of k n+ \ in P rest-frame we will use angles 1 <p i orienting t\ (also in P rest-frame). 
The Jacobian for this reparametrization of unit sphere equals 1 as wellj. 

Formula (TTST) can be iterated and provide a parametrization of the phase space with an ar- 
bitrary number of final state particles. In such a case, the question of orientation of the frames 
used to define the angles and the order of M,..„ integrations (consequently, the choice of limits 
for M,..„ integration), becomes particularly rich. Our choice is defined in ref. [2]. We will not 
elaborate on this point here, nothing new was introduced for the purpose of our study. 

If the invariant mass Mi... n is replaced with the photon energy defined in the P rest-frame, 
ky, then the phase space formula can be written as: 



dLips n+ i(P) 



kydkydcosQdty 



2(2ti) ; 



x dLips n (p—>ki...k„), (16) 



If we would have I photons accompanying n other particles, then the factor in square brackets 
is iterated. The statistical factor jj would complete the form of the phase space parametrization, 
similar to the formal expansion of the exponent. The last formula, supplemented with definition 
of frames with respect to which angles are defined is used to define the full kinematic config- 
uration of the event. From angles and energies (fcv ) of photons and also angles, energies and 
masses for other decay products, four-momenta of all final state particles can be constructed. 

If in formula (fT6l) instead of dLips n (p — > k\...k n ) one would use dLips n (P — > k\...k n ) the 
tangent space would be obtained. Then k n+ \ photon does not affect other particles' momenta 
at all, and thus has no boundaries on energy or direction. If this formula would be iterated 
then all such photons would be independent from one another as Energy and momentum 

constraints on the photon(s) are introduced with the relation between tangent and real n+l- 



4 The parametrizations of such a type, use properties of the Lorentz group in an explicit manner, in particular 
measure, representations and their products. That is why, they are useful, for event building Monte Carlo programs 
in phase-space constructions based on boosts and rotations. 

5 Let us point to another difference with respect to angles 9f , tyf used in formula (UJ (and also inrefs. [16], [21]). 
For example, if n = 2 then our dLips n= 2 phase-space is parametrized by the angles and (j) only. The two angles 
will define orientation of ky with respect to frame used for quantization of internal state of k\, this is different 
choice than the one used for 0^, even though = 9f . The Jacobian for the corresponding change of variables 
equals 1 also. 

6 Expression ( TTol l would be slightly more complicated if instead of photons a massive particle was to be added. 
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body phase-space. The formula defining one step in the iteration reads as followsB 

dUps n+1 (P^h. ..k n , k n+1 )= dLipsf tangent x C +1 , 

dLips+ l tangent = dkydcosQdty x dLips n (P -> h-h), 

{ki,...,k n +i} = T(ky,Q,$,{ki,...,k n }). (17) 

The W„ +l depends on details of T, and will be thus given later in formula (1231) . To justify (fTTT) . 
we have to convolute formula (fT5l) for Lips n+ \(P — > k\...k n ,k n+ \) with itself (for Lips n (p — > 

Lips n+ i(P -> ki...k n ,k n+l ) = L - n Lips 2 (P -> fc^+ip) x Lips n (p -> k\...k n ) 

in 

dM\ , , 

Lips n (p^k\...k n ) = ■ n Lips2(p^hp)xLips n -i(p -»fc 2 — fen) (18) 

271 

and use it also for Lips n (P — > ki...k n ): 

dM\ - _, 

Ups n {P -^k\...k n ) = - ■ n Lips2(P — > fcip ) x Lips n -\(p — > ki—K). (19) 

271 

Note that our tangent space of variables dkydcos Bd§ is unbounded from above and the limit 
is introduced by which is set to zero for the configuations outside the phase-sace. In 

principle, we should distinguish between variables like M2... n for invariant mass of &2 • • -k n and 
^2...w for invariant mass of k 2 ...k n , but in our choice for G n , G n+ \ below, M2... n = ^2...n and 
Mi „ is defined anyway for the n + 1-body phase space only. 

We direct the reader to refs. [1,2] for an alternative presentation. Let us remark that formula 
(flTl) is quite general, many options, motivated by the properties of the matrix elements, can be 
introduced. Generally the transformation T may differ from the choice to choice quite a lot. The 
most straightforward choice can be based on any n and n + 1 body phase-space parametrizations 
using invariant masses and angles (e.g. exactly as in TAUOLA [21] formulas 11 to 13). 

If 

G„ : Ml „,0i,4>i,Mf ... n ,e 2 ,4>2,...,6„_i,<k-i -> h---k„ (20) 

and 

Gn+l '■ ky, ■ kn i k n -\- 1 (21) 

then 

T = G n+l (k 1: Q^G-\k u ...,k n ))- (22) 

The ratio of the Jacobians (factors X 1 / 2 like in formula ([15]) . etc.) form the factor W„ +1 , which 
in our case is rather simple, 



w«+l = k _J_ x ^ /2 (W/ML„,ML„/MLJ 
y 2(2tt) 3 X^{\,m\/M\Ml JM 2 ) ' 



7 The {k\ , . . . , kn] can be identified with the event before the radiation of ky is introduced. 
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because of choice for G as explained in the Appendix. Note that ky = — 2m'"'" . There are 
additional benefits from such a choice. In all relations k.2 = Lk2, k n — Lk n and p' = Lp' 
common Lorentz transformation L is used. Transformation L is defined from k\,k\,p' ,p' and 
P; internal relations between four vectors k2...k n , (k2---k n ) are not needed. 
Formula (fT7l) can be realized algorithmically in the following way: 

1. For any point in n-body phase space (earlier generated event), described for example with 
the explicit configuration of four vectors ki...k n , coordinate variables can be calculated, 
using formula (|20l) . 



2. Photon variables can be generated according to Eq. (TT7T ). The weight W„ +l has to be also 
attributed. 

3. Variables obtained in this way from the old configuration and the one of a photon can be 
used to construct the new kinematical configuration for the n + 1-body final state. The 
phase- space weight, which is zero for configurations outside phase space boundaries, can 
be calculated at this point from (|17I23I) and finally combined with the matrix element. 

Here we have chosen two sub-groups of particles. The first one consisted of particle 1 alone, 
and the second, of particles 2 to n combined together. Obviously in the case of 2-body decays 
as discussed in this paper, there is not much choice when construction of the first photon is 
performed. 

By iteration, we can generalize formula (fTTT ) to the case of I photons and we write: 



1 ' 

dLips n+ i (P^k h ..kn, k n+ i.. .k n+ i) = — J\ 

i=l 



dk^d cos 9y,^Y,^+/- 1 



x dLips n (P — > ki ...k„ 



{ki,...,k n+ i} = T(ky p Qy p §y p T(...,T(ky 1 ,Qy 1 ,§y 1 ,{h,...,k n })...). (24) 

In this formula we can easily localize the tangent space for the multiple photon configuration. 
In this space, each photon is independent from other particles' momenta. Note that it is also pos- 
sible to fix upper boundary on ky j arbitrary high. Photons are independent one from another as 
well. Correlations appear later, thanks to iterated transformation T. The factors W^-_ l are cal- 
culated when constraints on each consecutive photon are introduced; the previously constructed 
ones are included in the n + i—1 system^. 

Of course, for the tangent space to be useful, the choice of the definition of T must be 
restricted at least by the condition {k\ , • • • k n } — > {k\, • ■ -k n } if all ky i — > 00 

It is important to realize that one has to choose matrix elements on the tangent space to 
complete the construction used in PHOTOS. The number and energies of photons will be gen- 
erated on the tangent space first. Regularization of (at least) soft singularity must be defined. 



Configurations of ky j which can not be resolved are replaced with the ones of that photon dropped out. 
9 In fact further constraints have to be fulfilled to enable presampling for the collinear singularities. Note that 
variables ky m ,Qy m ,$y m are used at a time of the m— th step of iteration only, and are not needed elsewhere in 
construction of the physical phase space; the same is true for invariants and angles M\ ,8i,<|>i, ... , 9„_i ,0„_ i — ► 
k\ . . . k n of (12012 U . which are also redefined at each step of the iteration. 
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Rejection, and event construction, is performed with the help of formula (TT71) for each consec- 
utive photon. It diminishes photon multiplicity with respect to the one defined for the tangent 
space. Of course, as rejection implements changes in phase-space density, a matrix element 
(with virtual corrections) of the physical space can be introduced as well. 

The treatment of the phase space presented here lies at the heart of the construction of 
PHOTOS kinematics, and was used since its beginning. It exhausts the case when there is only 
one charged particle in final state. For multiple charged particle final states new complication 
appear, because all collinear configurations need simultaneous attention, and not only the one 
along k\ direction. A presampler with multichannel generation is needed. In our case we follow 
the same method^ as explained in ref. [21]. 

In the standard version of PHOTOS, as published in [1,2], the following matrix element is 
used for single photon emission when there is only one charged particle in final state: 



M\ 2 photos = \A Bom \ 2 WT 3 M (25) 



where 



WTf d 



Ana 2(1-*) / m\ \ 1 + (3cos0 1 - \j 1 -w|cos6 



WT{WT 2 l + (\-x) 2 V l-(3 2 cos 2 0y 2 1-pcosG 



P 



mf 1 



M 2 (l-x) (l-x+(m 2 -m 2 )/M 2 ) 
M 2 n m 



ml = 4— — -L— — . (26) 



M M 2 -(m l +m 2 ) 21 K M 2 (l+m 2 /M 2 ) 



The old and lengthy approximation WT£ ld for WT3 implemented in standard PHOTOS is kept 
for compatibility with ref. [2]. The expression defining WT3 without approximations reads: 

\M\ 2 = \A Born \ 2 4Tza( qi ^-q— 
1 lexact 1 1 V h-ky Rky 



= |A Born | 2 4rox --y x 



Xfl,W]fl-?5)sin 2 e 



= |A Born | 2 4TOC ^ xWT 3 (P,k u k 2 ,k y ); (27) 

here, (k\ +k 2 ) 2 = x. In both, the standard and exact version of PHOTOS, the same phase- 
space parametrization and presampler for collinear and soft singularities are used. Together 
with WT 3 the following factors related to phase space WT\, WT 2 , contribute to the final weight 



Ml 2(f:) 2 (l + ^-y/2(i,^,^)cos0 
8 



10 We will omit details here, because for the two-body final states the complications manifest themselves only in 
the case of multiple photon generation, see Appendix. 



10 



implemented in routine PHOCOR of PHOTOS: 



WT l (P,k h k 2l ky) - J—; 



^jjj 2E y 



X ' X 



2(1 — cosGa/1 — mk) m 



The expression for WT{ can be deciphered from formula (fT71) and WT2 is related to presamplers 
for collinear and soft singularities. The factors WT\ and WT 2 are used in present paper in 
definition of WT£ ld (formula |26l) only. 

The combined effect of the virtual and real corrections on the total rate manifests through 

i-Total 

Ysss- The virtual corrections are included into PHOTOS through this factor. Let us point that 
the ratio of (l27l) and (1251) constitutes the basic element of upgrading PHOTOS functionality to 
the complete first order H The correcting weight can be chosen simply as: 

I of i 2 r Born 

wf= \M\exact (29) 
of 2 pTotal v ' 

\ Jn I PHOTOS 

For the standard version of PHOTOS the virtual corrections are required to be such that the total 
decay rate remains unchanged after complete QED corrections are included. 

In case of final states with two charged particles in PHOTOS the formula (|29l ) need to be 
modiffied with one of the following versions of the interference weight: 



wt = 



E 

1=1,2 







2 

exact 




2 

PHOTOS 



^Born 



pTotal 



WT, 



INT 



k\ .£ fc2-£ 

U/T , _ y qi h-h ^h.kj 

w 1 int — ~2 — : 



WT INT _ option 



Ji 



fci ,e fo.e 

yl^X C t 2 k 2 .k y 



1 



h WT x {PMM^WT 2 {PMM,ky) (30) 
The sum over two generation channels i = 1,2 related to emisson from q\ and q 2 is to be 
performed^. The form of WTj' NT results from the exact expressions, formulae (fh2~l) and d6]). 



11 When the option of multiple radiation is used in PHOTOS, the single photon emission kernel is iterated. This 
leads to some complications. 

12 Formulae d27T i and ( f25l l for \M \^ xact and |rW ^ PHO tos are f° r em i ss i° ns m case of single charge final state only, 
the interference weigt is to introduce the exact matrix element for process of two charged scalar final state. 
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However, phase space and multichannel presampler specific terms (|28l ) need to be discussed. 
Presence of J\ and J2 in interference weight is optional, but only for single photon radiation. 
The factor (J\ or J2) should cancel the WTi ■ WT2 term of the generation branch used for 
this particular event generation. In general, the absence of /1 2 terms is due to properties of the 
second order matrix element. For the time being analogies to the case of Z decay have to be 
used instead of the proof. 

Once we have completed the description of our internal correcting weight necessary for 
PHOTOS to work in the NLO regime, we will turn to the numerical results. 

4. Results of the tests 

The most attractive property of Monte Carlo is the possibility to implement selection criteria 
for the theoretical predictions that coincide with the experimental ones. Especially in the case 
of the final state bremsstrahlung presence of experimental cut-offs is essential, as they usually 
significantly increase the size of the QED effects. 

In this section we will concentrate however, on the following pseudo-observables, as used 
inref. [22,23]: 

• Photon energy in the decaying particle rest frame: this observable is sensitive mainly to 
the leading-log (i.e. collinear) non-infrared (i.e. not soft) component of the distributions. 

• Energy of the final-state charged particle: as the previous one, this observable is sensitive 
mainly to the leading-log (i.e. collinear) non-infrared (i.e. not soft) component of the 
distributions. 

• Angle of the photon with final-state charged particle: this observable is sensitive mainly 
to the non-collinear (i.e. non-leading-log) but soft (i.e. infrared) component of the distri- 
butions. 

• Acollinearity angle of the final-state scalars: this observable is sensitive mainly to the 
non-collinear (i.e. non-leading-log) and non-soft (i.e. non-infrared) component of the 
distributions. 

We will start our comparison for B -> n°K-(y) and PHOTOS running without improve- 
ments from the complete matrix element: the agreement looks good, see fig. [H and holds over 
the entire range of the distributions, even though densities vary by up to 8 orders of magnitude. 
Differences can hardly be seen. To visualize the differences, in fig. |2l the ratios of the distri- 
butions are plotted. Similar to what was seen in the tests for Z decays [14] local discrepancies 
may reach up to 15 % for cos0 aco //. > 0.5. Note however that those regions of the phase-space 
contribute at the level of 10 to the total decay rate. Once the matrix element is switched on, 
see fig. [31 where ratios of distribution are plotted, the agreement become excellent, even at a 

13 For example the form WTj^T-option is inappropriate for configurations when the first generated photon is hard 
and the second soft. 
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statistical level of 10 9 events. It was of no use to repeat the plots of the distributions with the 
corrected weight in PHOTOS, as the plots could not be distinguished from the ones of fig. \T\ 

Encouraged by the excellent performance in the case of the decay into final states with a 
single charged particle, let us now turn to decays into two charged mesons. To avoid accidental 
simplifications, we have selected final states with scalars of different masses 

Again, as can be seen from figs. @]and[5l agreement between PHOTOS using the stan- 
dard kernel and SANC is rather good, but some differences persist. Once the complete kernel 
is switched on, fig. |6l the agreement is quite amazing. In this case, the interference weight, 
and the multiple singularity structure of the pre- sampler Jacobians, formula (|30l) . were tested 
as well. Both versions WT INT and WT INT - option gave the same results for the case of single 
photon emission. However only the first version, WTj^j turned out to be consistent with ex- 
ponentiation. To complete the tests for multichannel emissions, final states with more than two 
massive decay products need to be studied, preferably for multi-photon radiation as well. 

Let us comment that not only the shapes of the distribution agree in an excellent manner 
between PHOTOS and SANC simulations, also the number of events with photons of energy 
below the certain threshold agreed better than 0.01 %, thus were consistent with each other 
within a statistical error of 10 9 event samples. The excellent agreement, presented in our paper, 
combined with other results published before, help to confirm that theoretical effects normally 
missing in PHOTOS are small, but if necessary can be introduced into the code. It is also 
important to note that the agreement provides powerful technical test of the generator. 

Finally, let us point out that early versions of the program, before 2004, were not reach- 
ing that level of technical sophistication. To establish it required a major effort. Kinemati- 
cal variables used in PHOTOS differ from those of SANC. The differences could arise due 
to technical problems, also if for example the Born-level events which are to be modified by 
PHOTOS would not fulfill energy-momentum conservation or particles momenta were not on 
mass-shell, at the numerical double precision level. This point must always be checked for 
every new installation of PHOTOS in an experimental environment. For that purpose we have 
collected numerical results, given in Table 2., for the cumulant of bremsstrahlung decay width: 
G(E test ) = r(£ fesf )/r Total , where T(E test ) denotes the decay width, integrated over energy car- 
ried by all bremsstrahlung photons combined up to maximum of E test . 

5. Summary 

This paper was devoted to the study of bremsstrahlung corrections in the decay of 5-mesons 
into pair of scalars of rather large masses. The results were presented in the analytical form and 
Monte Carlo simulations which were later compared. 

To quantify the size of the Next to Leading order effects normally missing in PHOTOS 
we have installed into the program the complete scalar-QED first order expression for the B 
decay matrix element. After modification, the differences between PHOTOS and the matrix 
element calculation embodied in SANC were below statistical error of 10 9 events for all of 
our benchmark distributions. Both PHOTOS and SANC were run at fixed first order without 
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exponentiation. The agreement provides a technical test of the simulations from both of the two 
programs as well. 

The improvement of the agreement due to the introduction of a correcting weight could 
come with a price. That was the case with the decay of Z. However because our 5-mesons are 
scalar, the complications did not materialize and a correcting weight can be installed to standard 
PHOTOS versions. On the other hand, introduced improvements are numerically small. Defi- 
ciencies of standard PHOTOS are localized in corners of bremsstrahlung phase space populated 
by photons of very high energies and angularly well separated from final state mesons. Those 
regions of the phase space weigh less than 0.005 to the total rate and differences in that region 
approach 20 % of their size, at most. The effects are thus significantly lower than 0.1 %, if 
quantified in units of the integrated B decay rate of a particular channel. Also, in those regions, 
the predictive power of scalar QED is rather doubtful. That is why we do not think it is urgent 
for users to change the PHOTOS correcting weight to enable the complete NLO, unless mea- 
sured form-factors become available. Contribution to the systematic error of PHOTOS due to 
incompleteness of the old kernel (with respect to scalar QED) does not depend on experimental 
cuts and is thus of no phenomenological importance for today. 

Our paper was not only focused on numerical results due to final state bremsstrahlung in B 
decays. Aspects of mathematical organization of the program for calculation of radiative cor- 
rections for B production and decay was discussed as well. Approximations used in PHOTOS 
affects matrix elements and no? phase space, which is treated exactly including all mass effects. 

Details of phase- space parametrization and other aspects necessary for implementation of 
NLO effects are collected for the first time. Generation of the phase space starts from the 
tangent space constructed from an eikonal approximation but used also for the hard photons, 
even of energies above the available maximum. In the second step, phase-space constraints are 
enforced. The method is similar to exclusive exponentiation [24]. 

Complete re-analysis of the final weight for decays into scalars was presented. Parts corre- 
sponding to matrix elements, phase space Jacobians and generator pre-samples were explicitly 
separated. Special care was devoted to mass terms. Analytic form of the single photon emis- 
sion kernel (i. e. matrix element with approximation) used in standard version of PHOTOS, 
was also explicitly given. That is why, the analysis presented here can be easily extended to 
other decay channels. It is the first time that we have presented such a study for particles other 

than elementary fermions and in the case where mass terms of order f -^r are not neglected. 
Our analytical calculations agree with the results of reference [19] exactly and could serve as a 
basis of our technical tests of the program. 

The numerical results collected in Table 2. can be used as a technical test of PHOTOS in- 
stallation in end-user simulation set-ups. We strongly recommend such tests to be performed. In 
these tests the agreement between PHOTOS and SANC (or simple semi-analytical expressions 
for higher order simulations) was significantly better than 0.1 % for all entries. 

In case of program operating for multiple photon radiation, energy momentum constraints 
are introduced for each consecutive photon, step by step, and conformal symmetry is not ex- 
ploited in that procedure. Deatails of phase-space parametrization used for multiple photon 
radiations were presented. In principle, for B decays and multiple photon radiation in PHO- 
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TOS, a similar level of agreement as in ref. [14] for Z decay is expected, but the appropriate 
reference distributions do not exist yet. In particular the second order, scalar-QED, matrix el- 
ement for B decays was not available for us. That is why we think, that the matrix-element 
related details of the program construction, necessary for implementation of NLL effects in 
general case, must remain delegated to the forthcoming work; probably to the times, when we 
will have at our disposal other second order matrix elements than just of Z — > //. At present, 
only dispersed (and for NLL based on analogies with Z decays only) results of refs. [2, 14, 25] 
are available for that purpose. 

PHOTOS used for decays of 5-mesons into scalars provides an example of program working 
for multiple emissions from both outgoing charged lines. It covers complete phase space, no 
special treatment is needed for the hard photon emission regions. Also mass terms have been 
included without any approximations. 

On the technical level it is worth mentioning that the NLO correcting weight of PHOTOS is 
used as an internal weight. All generated events remain weight 1, exactly as it was in the case 
of Z — > decay. 

In principle, if necessary, even complete higher order matrix elements (NNLO level) could 
be incorporated with the help of correcting weights. This interesting point definitely goes be- 
yond the scope of the present paper and also beyond the phenomenological interest for any 
foreseeable future. This is equally true for the possible extensions to simulations in QCD, 
which are also outside the scope of the paper. 
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Figure 1: Results from PHOTOS, standard version, and SANC for B — > k K~(y) decay are 
superimposed on the consecutive plots. Standard distributions, as defined in the text, and log- 
arithmic scales are used. The distributions from the two programs overlap almost completely. 
Samples of 10 9 events were used. The ultraviolet scale, /u uv , was chosen to leave total decay 
width unchanged by QED. 
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Figure 2: Results from PHOTOS, standard version, and S ANC for ratios of the B -> n ^ (y) 
distribution in figfflare presented. Differences between PHOTOS and SANC are small, but are 
clearly visible now. 
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Figure 3: Results from PHOTOS with the exact matrix element, and SANC for ratios of the 
B — > K°K~(y) distributions. Differences between PHOTOS and SANC are below statistical 
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Figure 4: Results from PHOTOS, standard version, and SANC for B — > n~K + (y) decay are 
superimposed on the consecutive plots. Standard distributions, as defined in the text and log- 
arithmic scales are used. The distributions from the two programs overlap almost completely. 
Samples of 10 9 events were used. The ultraviolet scale, /j uv , was chosen to leave total decay 
width unchanged by QED. 
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Table 2: Benchmark results for B decays into pair of scalar s: electromagnetic cumulative of 
decay width T(E test ) /r Total , where E test denotes the maximal energy which can be carried out by 
photons. The following input parameters were used: mg = 5279 MeV, m^o = 135 MeV, m % ± = 
139 MeV, m K o = 494 MeV, tn K ± = 498 MeV. Our results differ negligibly between standard 
PHOTOS and the one with exact matrix element that is why only one set of numerical results 
is provided. For each decay channel PHOTOS results of first, second and multiple photon 
radiation are to a good precision in the following proportion 1 — x : 1 — x + x 2 /2 : exp(-x), 
where x for each line of the Table is different; it depends on the decay channel and E test . To 
produce results for our table samples of 10 7 events were used. Statistical errors are thus at the 
level of the last significant digit for all the table entries. 



Appendix: Details and properties of the explicit phase space parametrization 

The formule (1 1712-4-1) are central in the definition of the PHOTOS algorithm and its phase- 
space parametrization. For the general description details were not important. In practice, they 
are nonetheless essential, and all angles masses (or energies) used in formulas (12012 II) must be 
specified. In particular all reference frames used in the definition of angles must be defined. 

We will start with the detailed description of parametrizations for two body and three body 
phase spaces, the latter one with an additional single photon which accompanies the final state 
of two massive (not necessarily of equal masses) objects. In both cases, the decay of an object 
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of the mass M and four momentum P is taken into account. The straigtforward extension for 
the parametrization of the multibody decay will be introduced with the help of the footnote, 
properties will be discussed later in the text. 

Our particular choice of the phase-space parametrization is of course motivated by the ne- 
cessity to regularize the infrared and collinear singularities. On the other hand, the definition 
itself does not need singularities to be exposed, or even to be present at all. To ease reading, let 
us point out that (at first) we will expect the collinear singularity to be present only when the 
photon becomes parallel to the direction of k\. Later, we will discuss the case when both final 
states (of momentum k\ and k^) are charged and thus the singularity may appear along the two 
directions. We will continue with the case when the photon accompanies a multi-particle/multi- 
charge final state and on the necessity to introduce several simultaneous parametrizations, to be 
used in Monte Carlo parallel generation channels, used at each step of iteration as defined in 
formula (f24|) . 

In the following eight points we define the angles used for the two-body phase-space parametriza- 
tion, and continue with the definition of phase-space variables of the two-body plus photon case. 

1. For the definition of coordinate system in the P-rest frame the x and y axes of the lab- 
oratory frame boosted to the rest frame of P can be used. The orthogonal right-handed 
system can be constructed with their help in a standard way. 

2. We choose polar angles 9i and §\ defining the orientation of the four momentum &2 in 
the rest frame of P. In that frame k\ and k\ are back to baci lli see fig. ©. 

3. The previous two points would complete the definition of the two-body phase space, if 
both k\ and k.2 had no measurable spin degrees of freedom visualizing themselves e.g. 
through correlations of the secondary decay products' momenta. Otherwise we need to 
know an additional angle §x to complete the set of Euler angles defining the relative 
orientation of the axes of the P rest-frame system with the coordinate system used in the 
rest-frame of k.2 (and possibly also of k\), see fig. ©. 

4. If both rest-frames of ki and £2 are of interest, their coordinate systems are oriented with 
respect to P with the help of 9i, §\, §x- We assume that the coordinate systems of k\ 
and £2 are connected by a boost along the k\ direction, and in fact share axes: z' Ti z" , 

5. Let us turn now to the three-body phase space parametrization. We take the photon energy 
kj in the rest frame of P, with its help we calculate: photon, k\ and &2 energies, all in 
k\ + k2 frame. 

6. We use the angles 0, (]), in the rest-frame of the k\ + k2 pair: angle is an angle between 
the photon and k\ direction (i.e. — z" )■ Angle (]) defines the photon azimuthal angle around 
z" , with respect to x" axis (of the &2 rest-frame), see fig. ©. 

14 In the case of phase space construction for multi-body decays ki should read as a state representing the sum 
of all decay products of P but k\ . 
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7. If all ki, &2 and k\ + k2 rest-frames exist, then the x-axes for the three frames are chosen 
to coincide. It is possible, because they are all connected by the boosts along the common 
z" direction, see fig. ©. The axes of k\ + &2 rest-fame are not drawn explicitly. 

8. To define orientation of &2 in P rest-frame coordinate system, and to complete construc- 
tion of the whole event, we will re-use Euler angles of k^'. §x, ©1 an d §i (see figs. [lOland 
[TT)) . defined again of course in the rest frame of P. 

That completes our definitions of parametrizations for the two-body and three-body phase 
spaces necessary to define transformation G of formulas (12012 II) . Before commenting on the 
properties of our parametrizations, let us note that these parametrizations were already used 
and defined in ref. [2], in all details, except for the function of (only implicitly introduced there) 
angle §x- F° r some readers, definition from that paper may be also easier to follow. 

Let us comment, now, on these properties of our parametrizations which are important for 
construction of the PHOTOS algorithm. 

a) The parametrizations of the two-body and three-body phase- space (photon included) are 
used for the explicit kinematical construction denoted by formula (fT7l) . We can replace the 
roles played by k\ and &2- This simple operation leads to a new phase-space parametriza- 
tion, which can be used in a second branch of the Monte Carlo generation. 

b) The phase-space Jacobians (factor W," +i of (flTT) ) are identical for the two branches; this 
factor is also never larger than ^Y 2 (27t) 3 ' 

c) Angle of the first branch coincides with 71 — of the second one. 

d) In the soft (ky — > 0) and collinear (0 — > or 71) limits, angles 0i , §\ , §x of the two branches 
converge to each other (in these limits they may differ by n or 2%). 

e) Properties (c) and (d) are convenient for our construction of the weights given by formula 
(|30l) , because they coincide with the similar properties of the exact matrix element. 

f) Thanks to property (b), also the first version of (1301) is exact. In fact, this first version is 
more suitable for multi-photon radiation, if first order matrix element is used only. This 
observation required comparisons with second order matrix elements [14]. The choice of 
k.2 (or ki) direction to define 0i, §\, rather than k\, was also motivated by the properties 
of the decay matrix elements. 

Let us present now some further observation^, which go beyond NLO corrections for the 
processes discussed in ref. [14] and in the present paper, and point to applications for a multi- 
body /multi-charge final states. 

15 Note, that the approximations to be discussed in the following points, result from matching kinematical 
branches and affect the way how phase-space Jacobians are used in ( f3Qb . The full phase space remain covered, as 
is the case of formule (117124b denoting exact phase-space parametrization of single and multi-photon final state. 
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• Property (d) extends to more than two-body decays, and also to cases when there are more 
than two charged particles present in the final state. The relation between angles 0i, §i, 
§X °f me distinct branches is more complex, but in the discussed limits still independent 
from and (|). 

• Extended property (d) and property (e) enable the use of (1301) for multi-photon radiation; 
this also holds in the case when more than two charged particles are present in the final 
state. 

• That is why, in the case of two-body decays (plus bremsstrahlung photons), such type of 
phase-space treatment is sufficient for the NLO precision. 

• For the NNLO precision, in matching of the two mappings for the collinear singularities^! 
another factor of the type X l ^ 2 {...) /?i 1//2 (...) would have to be included in W£ +l of for- 
mulas (1 171301) . In fact in such a case the exact multi -photon phase space parametrization 
would be preserved. 

• For each additional charged decay product present in the final state, still another factor of 
the type l 1/2 (...) A 1/2 (-..) is needed in W£ +l to assure multichannel generation with the 
exact treatment of the phase space. 

• Even without future refinements (as explained in the previous two points) our phase space 
parametrization is sufficient for NLO and NLL precision for the two-body (two-charges) 
decays, accompanied with arbitrary number of photons. In a general case, when more 
than two charged particles are present in final state, such phase space parametrization 
remains sufficient for LL only, even though also in this case, the full multi-photon phase 
space is covered. At present, the resulting precision is sufficient and does not justify those 
easy changes. 

• In our choice of phase space parametrization (point 1), we have dropped some details, the 
choice of x, y, z axes of the P rest-frame were not specified. Indeed, for the decay of a 
scalar object, such as that discussed in the present paper, every choice is equivalent. In 
general, it is not the case. Already in case of the Z boson decay, the choice of the z axis 
parallel to the direction of the incoming beam of the same charge as &2 is advantageous, 
see ref. [14] where the process e + e~ — > Z — > l + l~n(y) was studied. In this case the 
direction of the incoming beam coincides with the spin state of Z, and the choice simplify 
expression for matrix element. 

• Finally, we may consider extending our method beyond the decays. In the case of the 
?-channel processes, or initial state radiation, particular choices for the frames and angles 
would be also essential for the construction to match the structure of the matrix element 
singularities. Some hints into that direction can be seen already now, refs. [25,26]. 



Such matching is necessary for the two branches of the generation, used to presample collinear singularities 
along the directions of k\ and k%, to be used simultaneously in construction of each event. 
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Figure 7: The angles 9i, <|)i defined in the rest-frame of P and used in parametrization of two- 
body phase-space. 




Figure 8: Angle §x is also defined in the rest-frame of P as an angle between (oriented) planes 
spanned on: (i) k\ and z-axis of the P rest-frame system, and (ii) k\ and /'-axis of the k.2 rest 
frame. It completes definition of the phase-space variables if internal orientation of k\ system 
is of interest. In fact, Euler angle §x is inherited from unspecified in details, parametrization of 
phase space used to describe possible future decay of £2 (or k\). 




Figure 9: The angles 0, § are used to construct the four-momentum of ky in the rest-frame of 
k\ +&2 pair (itself not yet oriented with respect to P rest-frame). To calculate energies of k\, ki 
and photon, it is enough to know m\, m,2, M and photon energy ky of the P rest- frame. 
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Figure 10: Use of angle ty x in defining orientation of k\, &2 and photon in the rest-frame of P. 
At this step only the plane spanned on P frame axis z and &2 is oriented with respect to k2 x x" 
plane. 




Figure 11: Final step in event construction. Angles 0i, §\ are used. The final orientation of Ic2 
coincide with this of ^2- 
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